EVOLUTION OF FERMIONIC SYSTEMS AS 
AN EXPECTATION OVER POISSON PROCESSES 
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C^") ' We derive an exact probabilistic representation for the evolution of a Hubbard 

C _ ) model with site- and spin-dependent hopping coefficients and site-dependent in- 

teractions in terms of an associated stochastic dynamics of a collection of Poisson 
processes. 

\ In a recent paperB we provided an exact probabilistic expression for the real 

^ ■ time or the imaginary time evolution of a Fermi system, in particular a Hubbard 

model, in terms of an associated stochastic dynamics of a collection of Poisson 
processes. Here, we generalize the resulti!l to the case of a Hubbard model with 
site- and spin-dependent hopping coefficients and site-dependent interactions. This 
situation is of interest if disorder is presenta. 
Let us consider the Hubbard Hamiltonian 



|A| |A| _^ |A| 
i=l j=i £T=T4. <=1 



where Ac Z d is a finite <i-dimensional lattice with cardinality |A|,{1,...,|A|} some 
total ordering of the lattice points, and Ci a the usual anticommuting destruction 
operators at site i and spin index a. Note that the Hamiltonian ([!]) also allows 
for spin-dependent site energies —2rju a . We are interested in evaluating the ma- 
trix elements (n'\e~ lHt \n) where n = (ni^mi, ■ ■ ■ ) n \A\1> n \A\l) are the occupation 
numbers taking the values or 1. Since the total number of fermions per spin 
component is a conserved quantity, we consider only configurations n and n' such 

that J2[=i n 'ia — Sl=i n iv f° r CT = TI- I 11 the following we shall use the modulus 2 
addition n © n! = (n + n') mod 2. 

Let = < i < j < |A| : r]ij a ^ 0} and |r ff | its cardinality. By 

introducing 

X ija (n) = (n 8 l f(T © l j<7 \4 a c ja + c] a c ta \n) 

= \n 3a {n la © 1) - n la (n 3a © 1)] , (2) 
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where l ia = (0, . . . , 0, l ia , 0, . . . , 0), and 

|A| |A| 

V(n) = (n|i?|n) = E^^t^l - E E 2-q UlJ n la , (3) 

i=l <=1 ct=U 

the following representation holds 

(n'|e-^ 4 |n) =E(<5 n ^„ eNt X*) (4) 

M* = ex P { E E / lo S \inu,p£*ij*(n © N s )] diV?- CT 
L -=Ti(M)er.- / I°.*) 

-i [ V(n(BN s )ds+ E ^4" (5) 

Here, {iV*^}, (i, j) e and a =tl, is a family of [T-f- ] + |TJ independent Pois- 
son processes with parameters pij„ and N* = (NlpNf^, . . . , N^ A ^, N* A ^) are 2|A| 
stochastic processes defined as 

NL= E N h*- ( 6 ) 

(l,i)Gra-' i—k or j — 

We remind that a Poisson process AT* with parameter p is a jump process char- 
acterized by the probabilities P (N t+S — iV* = k) = (ps) e~ ps /k\. Its trajectories 
are piecewise-constant increasing integer- valued functions continuous from the left. 
The stochastic integral J dN* is just an ordinary Stieltjes integral 



/ f(s,N s )dN s = E f( s k,N Sk ), 

J\o.t) ,.. : ^ 



'[°.*) k: s k <t 

where are random jump times having probability density p(s) — pe~ ps . Finally, 
the symbol E(. . .) is the expectation of the stochastic functional within braces. 
The matrix elements (n'\e~ lHt |n) obey the ODE system 

|<n'|e-"» = -iE( n 'l^l n ")(n"|e-^|n) (7) 

n" 

with initial condition (n' | e — tff * | n) | t _ Q = <5 n 'n- One may check that (fl||5|) is indeed 
solution of (^) by applying the rules of stochastic differentiation. We have 

E (^n' : nffiN'+ t!t - A ^' +dt ) 

= E E (ll II 5n',n»e d N«e / [«+^ log [ i "«''' , «^*^ (nffiNS) ] diV ^ 
n" ^cr=Tl (i,i)er„ 

x e -iV(neN*)rft+E CT=n E (i;3) er. PH.* S a „ fBmK *M*) (8) 



For the Markov property, the expectation of the factors containing the stochastic 
integrals in the interval [t, t + dt] can be taken separately. By expanding each one 
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of them over all possible numbers of jumps of the Poisson processes as 

= <W e° e-^^' + ^^ei^ei,, e *-^ A -( n ® Nt )] < ",,„•*'" • ... 

= <5n',n" + [^n"ffil ta ffil ie ,%iaPyffAjjv(n©N*) -<^n',n"] /OijcrCft + O (dt 2 ) , 

to order dt we obtain 

*n',n» + 2J E 5 n',n"©li<,ei,V%iCTAy CT (n")dt 

-£ n /, n »iV(n")eft 

Finally, we rewrite this relation as 

(V.neN'M') = E (5 n ,, neNt+dt A^ t+dt ) - E (* n /, nfflN .A<*) 

= -i^2(n'\H |n")E (V^n^*) A. (10) 
n" 

It is clear that the fcrmionic nature of plays no special role in the above derivation 
which holds for any system described by a finite Hamiltonian matrix. 

In order to construct an efficient algorithm for evaluating we start by 

observing that the functions Aij <T (n © N s ) vanish when the occupation numbers 
rii a © Nf a and Uj a © N* are equal. We say that for a given value of a the link 
ij is active at time s if A-y CT (n © N s ) 7^ 0. We shall see in a moment that only 
active links are relevant. Let us consider how the stochastic integral in (^J) builds 
up along a trajectory defined by considering the time ordered succession of jumps 
in the family {N^ a }. The contribution to the stochastic integral in the exponent 
of (||) at the time of the first jump, for definiteness suppose that the link iijx with 
spin component <Ti jumps first at time si, is 

log [iVi 1 jicxPl 1 j 1 ai^*ihn(p © N S1 )] °( t ~ s i)> 
where N Sl = due the assumed left continuity. Therefore, we obtain a real finite or 
—00 contribution to the stochastic integral depending if the link was active 

or not at time 0. If s\ > t we have no contribution to the stochastic integral from 
this trajectory. If si < t a second jump of a link, suppose with spin component 
<72, can take place at time S2 > s\. The analysis can be repeated by considering an 
arbitrary number of jumps. Of course, when the real part of the stochastic integral 
is —00, i.e. when some A = 0, there is no contribution to the expectation. The other 
integral in (JsJ) is an ordinary integral of a piecewise constant bounded function. 

From the above remarks it is clear that the only trajectories to be considered 
are those associated to the jumps of active links. We start by determining the set 
of active links A\ C r ff , a =TI, in the initial configuration n assigned at time 0. 
This is done by inspecting the occupation numbers of the sites according to the 
rule that the link ij is active for the spin component a if Uj CT + rij a = 1. Then, we 
extract the jump time Si according to the probability density 

p( s )=E E ^e- E "^<-^i ft - s . (11) 
<r=U (i,j)eAi 
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E (<5n',n®N*+ dt -M t+dt ) = E 



E(<S n 



PN 



:M* 



(9) 



Finally, we choose the jumping link in the set A\ with probability 

p.. - (12) 

The contribution to M 1 at the time of the first jump is therefore, up to the last 
factor which appears in (|5|), 

+ e -iV(n©N f )t e -E< T=Ti E( s , 3)e r„\Ai P^-* 0( Sl _^ 

where exp ^— X)o-=tJ, S(i j)er VA 1 Pijc^j IS the probability that the non active links 
do not jump in the time interval s. The contribution of a given trajectory is obtained 
by multiplying the factors corresponding to the different jumps until the last jump 
takes place later than t. For a given trajectory we thus have 

m* = nKi^^r4^ A ^^( n ® Nsfc ) e[a " iv(neN8 ' e)i(s ' ! " s ' ! - i) 

k>l 

+e [C»-iV(nffiN*)](t-a»_0 0(sfc _ t )j _ (13) 

Here, Cfe — So-=|x E(j j)eA fc Ptf" wri ere is the set of active links with spin a 
in the interval (s/-_i,Sfe] and So = 0. Note that the last exponentially increasing 
factor in (|[) cancels out in the final expression of M*. The analogous expression 
of Mr for imaginary times is simply obtained by replacing 7/ — > — ir\ and 7 — > —17. 

In principle, the algorithms parametrized by different are all equivalent 
as (^-||) holds for any choice of the Poisson rates. However, since we estimate 
numerically the expectation values with a finite number of trajectories, this may 
introduce a systematic error. It can be shown that the best performance is obtained 
for the natural choice pij„ ~ r\ij a independently of the interaction strength 7. 
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